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UNDERWATER  TARGET  TRACKING  USING  MEASUREMENTS  FROM  A  LINEAR  ARRAY 
OF  OMNIDIRECTIONAL  ACOUSTIC  ELEMENTS 


1.  INTRODUCTION 

The  tracking  problem  discussed  herein  is  defined  by  a  moving  observer 
in  an  underwater  environment  which  uses  measurements  from  a  trailing 
linear  acoustic  array  to  estimate  the  initial  position  and  velocity  of  a 
target.  The  assumptions  made  in  discussing  this  problem  are  as  follows. 
The  array  is  a  multispot  linear  array  of  uniformly-spaced  omnidirectional 
elements  operating  in  a  passive  mode;  it  is  considered  to  be  an  inflex¬ 
ible  body  which  is  rigidly  connected  to  the  observer's  stern.  The 
target-observer  kinematics  are  described  in  a  three-dimensional  Cartesian 
coordinate  system.  Both  the  velocity  and  depth  of  the  target  are 
constant  throughout  the  problem.  The  observer  maneuvers  to  improve  the 
quality  of  the  estimate  and  these  maneuvers  are  generated  by  impulses  in 
acceleration  and  are  restricted  to  a  horizontal  plane  at  a  depth  that  may 
differ  from  the  target's.  The  underwater  environment  is  modeled  as  an 
infinite,  homogeneous  medium  with  the  target  being  the  only  noise  source 
in  this  medium.  Finally,  it  is  assumed  that  the  linear  array  is 
operating  in  the  far-field  region,  so  that  only  plane-wave  prop¬ 
agation  is  considered. 

Various  versions  of  this  target  tracking  problem  are  produced  by 
changing  the  type  of  sensor,  altering  the  restrictions  on  the  environ¬ 
ment,  and  relaxing  the  restrictions  on  target  motion  and  other  assump¬ 
tions.  One  well-studied  version  is  the  bearings-only  case,  where  both 
the  tai ^et  and  observer  are  operating  at  the  same  depth  and  a  bow-mounted 
sonar  measures  the  azimuthal  angle  between  north  and  the  line-of-sight  to 
the  target.  In  other  words,  by  steering  in  a  horizontal  plane,  this  sonar 
isolates  a  line  segment  on  which  the  target  lies.  In  contrast  with  this, 
the  linear  array  measures  a  conical  angle  between  the  array’s  longi¬ 
tudinal  center  line  and  the  surface  of  a  cone  on  which  the  target  lies. 


The  measurement  equation  describing  this  conical  angle  is  a  nonlinear 
function  of  the  target  parameters:  velocity  and  initial  position. 
Consequently,  an  iterative  least-squares  estimator  is  used.  Although  the 
extended  Kalman  filter  is  also  applicable,  it  was  not  tried  because  of 
its  erratic  behavior  with  the  Cartesian  model  of  bearings-only  target 
motion  analysis,^1  which  has  many  of  the  same  properties  as  the  omni¬ 
directional  linear  array  problem.  Also,  the  additional  experience  with 
iterative  techniques  may  prove  useful  in  treating  other  problems  where 
past  measurements  from  secondary  sensors  must  be  incorporated  into  the 
current  estimate.  The  measurement  noise  is  assumed  to  be  Gaussian  even 
though  the  iterative  least-squares  estimator  does  not  require  any 
information  about  the  noise  structure. 

2 

The  modified  Gauss-Newton  method  is  the  particular  iterative 
technique  used  in  this  study.  However,  its  unmodified  version,  the 
Gauss-Newton  method,  is  also  discussed  and  some  results  obtained  with  it 
are  cited.  The  performance  index,  P I 1 ,  for  this  nonlinear  least-squares 
problem  is  the  sum  of  the  squared  measurement  residuals,  which  is  mini¬ 
mized  by  iteratively  estimating  the  target  parameters.  The  previous 
estimate  defines  a  nominal  target  track  which  is  used  to  formulate  a 
linearized  least-squares  problem  that  is  solved  by  a  Householder  trans¬ 
formation.-^  This  solution  updates  the  previous  estimate,  and  the 
procedure  repeats  until  the  terminating  criteria  are  satisfied.  The 
modified  Gauss-Newton  method  is  used  because  of  its  elementary  structure, 
and  the  Householder  transformation  is  used  because  of  its  numerical 
accuracy. 

The  next  section  of  this  report  develops  the  Cartesian  signal  model 

required  by  the  estimation  algorithm.  This  is  followed  by  a  section  that 

describes  the  Gauss-Newton  method,  the  modified  Gauss-Newton  method,  and 

4 

the  Householder  transformation.  The  Levenberg-Marquardt  method  is 
also  mentioned  because  of  its  success  with  singular  or  ill-conditioned 


2 


problems.  The  experimental  results  are  given  in  the  fourth  section,  and 
the  conclusions  and  recommendations  for  further  work  are  given  in  the 
last  section. 

Appendix  A  shows  the  Jocobian  matrix  required  in  the  Gauss-Newton 
method,  and  appendix  B  shows  a  scheme  for  selecting  the  step  size  in  the 
modified  Gauss-Newton  method.  The  number  of  computer  operations 
(additions,  multiplications,  divisions,  square  roots,  and  trigonometric 
functions)  required  for  one  iteration  of  the  modified  Gauss-Newton  method 
is  given  in  appendix  C.  A  listing  of  the  computer  code  used  to  define 
the  experiments  is  contained  in  appendix  0. 
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2.  SIGNAL  MODEL 


The  signal  model  is  defined  by  the  measurement  process  and  the 
target-observer  kinematics,  and  both  of  these  are  described  in  a  three- 
dimensional  Cartesian  coordinate  system  in  this  section.  The  complexity 
of  this  description  is  reduced  because  neither  the  target  nor  the 
ooserver  changes  depth.  A  further  simplification  results  from  the 
assumption  on  the  observer's  maneuvers;  i.e.,  instantaneous  changes  in 
speed  and  heading  are  permitted.  Consequently,  the  center  lines  of  the 
observer  and  the  linear  array  are  always  coincident.  This  gives  a 
shorter  convergence  time  when  compared  with  the  results  from  a  signal 
model  that  uses  a  realistic  description  of  the  maneuvers. 

Figure  1  shows  the  two-element  model  of  the  linear  array  used  in  this 
report.  The  distance  between  elements  is  b,  and  the  distance  traveled  by 
the  plane  wave  as  it  moves  from  the  first  to  the  second  element  is  d. 

The  direction  to  the  target  is  given  by  the  angle  a,  which  should  be 
defined  at  the  midpoint  between  the  two  elements.  This  error  is 
neglected,  as  is  the  offset  between  the  positions  of  the  array  and  the 
observer. 

Inspection  of  figure  1  shows  that 

cos  a  =  d/b  *  c  t^/b,  (2-1) 

where  c  is  the  nominal  speed  of  sound  in  sea  water  and  t^  is  the  wave 
front  time  delay  between  elements.  The  observer  reconstructs  the  plane 
wave  by  combining  the  electrical  signals  from  these  two  elements  after 
first  inserting  the  time  delay,  t^,  into  the  signal  path  of  the  second 
element.  Further  inspection  of  figure  1  shows  that  a  and,  consequently, 
t^  are  unchanged  when  the  direction  of  propagation  line  is  rotated 


about  the  array's  longitudinal  center  line.  Hence,  a  time  delay  measure¬ 
ment  defines  the  surface  of  a  cone  on  which  the  target  lies.  This  is 
shown  in  figure  2. 

In  this  report  t^  rather  than  a  is  defined  as  the  measured  variable 
because  cos  a  is  the  required  input  signal  for  the  estimation  algorithm 
presented  in  the  next  section.  Assuming  that  a  is  the  measured  variable, 
the  noise  level  on  t^  is  established  from  equation  (2-1).  The  noisy 
version  of  this  equation  is 

cos(a  ♦  va)  ,  c(tj  ♦  vt)/b  (2-2) 

where  va  and  vt  are  the  noise  terms.  Expanding  the  cosine  term  and 
simplifying  the  result  gives 

vt  *  va  Wc)  sin  a  (2-3) 

2 

Assuming  that  v  *  N(0,  <j  )  and  that  at  any  point  in  time  a  is 
a  a 

unknown  but  constant  over  an  ensemble  of  experiments,  then  vt  - 
N(0,  a^)  where 

ot  =»  b  |  sin  a|aa/c.  (2-4) 

Hereafter, 

®t  *  b  (2-5) 

is  used.  This  approximation  leads  to  a  pessimistic  result. 

The  noisy  time  delay  is  defined  as  tm;  i.e., 

»  h  *  v  (2-6) 
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Figure  1.  A  Plane  Wave  Arriving  at  a  Two-Element 
Omnidirectional  Linear  Array 


Figure  2.  The  Surfaces  of  Four  Cones,  Each  for  a  Different  Time  Delay 


where  the  t  subscript  on  v  has  been  dropped  for  brevity.  Solving 
equation  (2-1)  for  t^  and  substituting  this  result  into  equation  (2-6) 
yields 

tm(k)  =  (b  cos  o(k))/c  +  v (k )  (2-7) 

where  k  denotes  the  time  t  =  kh  at  which  the  measurement  is  taken  with  a 
fixed  sampling  period  h  seconds.  Equation  (2-7)  defines  the  measured 
time  delay  as  a  function  of  the  conical  angle  a.  To  define  it  as  a 
function  of  the  target's  position  variables,  a  Cartesian  coordinate 
system  is  introduced. 

A  (R^,  R-j-i  Rz)  coordinate  system  fixed  at  the  array  with 
along  the  center  line,  Rz  in  the  depth  direction,  and  orthogonal 
to  both  Rj^  and  Rz  is  shown  in  figure  3.  Inspection  of  this  figure 
shows  that 

cos  a  =»  R^/R^  (2-8) 

where 

R$  ,  (Ru2  ♦  RT2  +  Rz2)1/2  (2-9) 

is  the  slant  range  to  the  target.  Substituting  equation  (2-8)  into 
equation  (2-7)  gives 

tm  (k)  =  b  R^  (k)/c  R$(k)  +  v(k),  (2-10) 

which  shows  the  measurement  as  a  function  of  the  target's  relative 
position  variables.  Also,  from  figure  3,  it  is  easily  shown  that 


R^2  -  R$2  COS2  a  »  0 


(2-11) 


defines  the  surfaces  of  the  cone  on  which  the  target  lies. 

Because  the  (R^,  Ry)  plane  rotates  when  observer  changes 
course,  these  variables  must  be  transformed  to  a  north-referenced 
coordinate  system  where  the  constant  velocity  target  is  defined.  The  new 
variables  are  Ry  (which  points  north)  and  Rx  (which  is  orthogonal  to 
both  Ry  and  R£).  The  (Rx,  Ry)  and  (Ru,  Ry)  planes  are 
coincident  but  differ  by  a  rotation,  and  Rz  is  the  same  in  both  coor¬ 
dinate  systems.  Both  of  these  coordinate  systems  are  fixed  at  the  array. 

This  rotation  is  defined  by  the  transformation  shown  in  the  next 
equation. 


sin  cn  Rj 


Ry  -sin  cn 


cos  cn  Rn, 


(2-12) 


where  c^y  is  the  angle  betwen  the  Ry  and  the  R^  axes,  measured 
positive  clockwise  from  the  north  reference.  Because  this  is  an  orth¬ 
ogonal  transformation,  Rs  is  also  given  by 

R$  ,  (Rx2  ♦  Ry2  ♦  Rz2)1/2  .  (2' 


(2-13) 


Figure  4  shows  these  two  Cartesian  coordinate  systems,  the  array,  and  the 
observer's  heading  angle. 

This  transformation  Is  expressed  in  terms  of  the  observer's  heading 
angle,  ip,  by  substituting 


Cjj  »  \p  -  180* 


(2-14) 


T 
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Figure  3.  The  (Rn,  Rf,  Rz)  Coordinates  System  Fixed  at  the 


into  equation  (2-12).  After  first  multiplying  equation  (2-12)  by  the 
inverse  transformation,  this  gives 


s'  V 

s 

f  s. 

Rt 

■ 

-cos  4/ 

sin  ^ 

Rx 

R11 

—sin 

-cos  ^ 

Ry 

✓ 

Using  the  last  equation  in  the  matrix  equation  (2-15)  to  eliminate 
from  equation  (2-10)  gives 

tm(k)  «=  -b(Rx(k)sin  $  (k)  +  Ry(k)cos  i>  (k))/c  R$(k)  +  v(k)  (2-16) 

where  R$  is  given  by  equation  (2-13).  Equation  (2-16)  shows  the 
measured  variable  as  a  function  of  relative  position  variables  which  are 
defined  in  a  north-referenced  system. 

Inspection  of  equation  (2-16)  shows  that  the  sign  on  Rz  cannot  be 
resolved  by  this  measurement.  This  was  demonstrated  by  simulation, 
although  this  elementary  experiment  is  not  included  in  this  report. 


For  a  constant  velocity  target,  these  relative  position  variables  are 
defined  by 


Rx00 

RY(k) 

Rz(k) 


xtkh  +  Rxt(°)  -  RX0(M 

( 2—17  a ) 

ytkh  -  Ryt(0)  -  Ry0(k) 

(2— 17b ) 

izt(0) 

( 2—17  c ) 

where  Vxt,  Vyt,  Rxt(0),  Ryt(0),  and  Rzt(0)  define  the  target 
paramters  in  another  north-referenced  coordinate  system  which  is  fixed  in 


space  and  which  at  t  *  0  is  coincident  with  the  north-referenced  system 
fixed  to  the  observer.  Also,  RxQ(k)  and  Ry0(k)  define  the  observer's 
position  in  this  stationary  coordinate  system. 

Equation  (2-16)  defines  the  measurement  process  and  equation  (2-17) 
defines  the  target-observer  kinematics.  Together,  these  equations  define 
the  signal  model  required  by  the  estimation  algorithm  presented  in  the 
next  section. 


3.  ITERATIVE  LEAST-SQUARES  TECHNIQUE 


Because  the  measurement  equation  (2-16)  is  nonlinear,  an  interative 
least-squares  algorithm  which  batch  processes  a  set  of  measurements  is 
selected  to  estimate  the  target  parameters.  The  extended  Kalman  filter 
could  be  used  for  this  problem,  but  it  was  not  tried  because  its  perfor¬ 
mance  with  a  Cartesian  signal  model  of  the  bearings-only  problem  is 
suspect . 

In  this  least-squares  problem  the  modified  Gauss-Newton  method  is 
used  to  iteratively  minimize  the  sum  of  the  squared  measurement 
residuals.  At  each  iteration  the  Gauss-Newton  equations  are  solved  by  a 
Householder  transformation,  and  this  yields  an  increment  or  step  which  is 
added  to  the  previous  estimate.  The  modified  Gauss-Newton  method  first 
multiplies  this  step  by  a  positive  scalar,  which  is  selected  at  each 
iteration.  The  Levenburg-Marquardt  method  is  also  mentioned  in  this 
section;  however,  it  did  not  demonstrate  any  advantage  in  the  experi¬ 
ments,  so  its  use  is  not  recommended  at  this  time.  These  experiments  are 
not  contained  in  this  report. 

The  sum  of  the  squared  measurement  residuals  which  defines  the 
performance  index  PI1  is  given  by 


PIl(x) 


(f(x)| 


fT(x)f(x) 


K-l 

z 

k  *0 


fk2(x) 


where 


K  >  5,  xT  =  (Vxt,  V^,  Rxt  (0),  R^O),  Rzt(0)), 
fT(x)  -  (fQ(x),  fj(x),  ....  f<_1(x)),  and 


(3-1) 


V*)  -  (Rx(k )sin  (k)  +  Ry(k)cos  <l>  (k))  ♦  c  R$(k)tm(k)/b  (3-2) 
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In  this  section  x  represents  the  estimate  of  the  target  parameters. 
Equation  (3-2)  defines  the  kth  measurement  residual.  With  K  finite,  it 
is  expected  that  the  true  target  parameters  can  minimize  (3-1)  only  when 
the  measurements  are  noiseless. 

A  general  iterative  procedure  for  minimizing  PI1  is  given  by 

xi+1  -  xi  +  9i  (3-3) 

where  x^  is  the  (i+l)th  estimate,  x^  is  the  ith  estimate,  and  g^ 
is  the  ith  step.  In  theory,  the  initial  estimate  xQ  is  usually 
selected  as  zero.  Unfortunately,  this  did  not  work  well,  so  the  first 
measurement  is  used  to  aid  in  the  selection  of  Xq.  The  procedure 
terminates  in  a  finite  number  of  iterations  when  no  further  significant 
improvement  in  PI1  occurs  and/or  when  I  U1+i  -  x^  ||  is  less  than 

max 

some  threshold.  A  requirement  for  convergence  of  this  procedure  using 
the  modified  Gauss-Newton  method  is  given  in  reference  5.  Although 
convergence  proofs  are  beneficial,  they  do  not  guarantee  the  performance 
of  the  algorithm.  Even  though  convergence  is  assured,  double  precision 
arithmetic  and/or  an  excessive  number  of  iterations  may  be  needed. 

The  linearized  least-squares  problem  which  is  solved  at  each 
iteration  is  derived  from  the  following  Taylor  series  expansion  of  PI 1 : 

PIl(xi+1)  =.  PIl(xi+g.)  =  PIl(xi)  +  ’?pnT(xi)gi  +  higher-order  terms  (3-4) 

where  v  is  the  gradient  operator.  If  xi+^  yields  a  local  minimum,  then 
vPU(x.+1)  =  0,  and  the  gradient  of  equation  (3-4)  gives 


where  the  higher-order  terms  have  been  neglected.  The  step  g^  to  the 
next  estimate  is  given  by  the  solution  to  the  matrix  equation 


72PIl(xi)gi  »  7pIl(Xi) 


(3-6) 


where 


vPIl(x)  -  9  fT(x)  f(x)/9X  >  2  JT(x)  f(x) 


(3-7) 


and 


.  T  K-l  , 

vZPIl(x)  -  2  JT(x)  J(x)  +  2  J  (/  f  (*))  f.  (*).  (3-8) 

k»0  K  * 

J(x),  the  K  by  5  Jacobian  matrix  of  f(x),  is  given  by 

J(x)  .  a  f(x)/ax  »  [3  fk/ a  x-j],  (3-9) 

which  is  given  in  appendix  A.  This  appendix  shows  that  the  last  column 
of  J(x)  is  zero  if  the  estimate  of  the  target’s  depth  is  zero.  This 
singularity  did  not  present  any  difficulties  for  the  experiments  reported 
here;  it  will  be  discussed  in  detail  in  a  future  report. 

The  computer  burden  is  reduced  if  the  last  term  in  equation  (3-8)  is 
replaced  by  an  approximation  which  does  not  require  any  second  deriva¬ 
tives.  Because  this  term  depends  upon  the  residuals,  the  approximation 
usually  introduces  negligible  error  near  the  solution.  Various  iterative 
techniques  result  from  different  approximations.  Using  equation  (3-8)  as 
stated  gives  Newton's  method,  while  dropping  the  last  term  gives  the 
Gauss-Newton  method.  In  addition,  setting  ?2PIl(x)  »  I,  an  identity 
matrix,  gives  a  steepest  decent  step  g1  „  _7PII(xi).  The 
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Levenberg-Marquardt  method  provides  a  mix  of  these  latter  two  techniques. 

Matrix-updating  methods,  such  as  the  Davidon-Flectcher-Powel 1  method, 

2 

approximate  the  Hessian  matrix,  7  PI,  by  a  successive  addition  of 
low-rank  matrices. 

The  convergence  rate  for  Newton's  method  is  second  order,  while 
reference  5  states  that  at  most  the  convergence  rate  for  the  Gauss-Newton 
method  is  linear  unless  P I 1 ( x  .  )  =  0,  in  which  case  it  is  ultimately 
second  order. 

Dropping  the  last  term  in  equation  (3-8)  gives 

JT(xi)  J(Xi)  gi  =  -JT(x.)  f(xi),  (3-10) 

which  defines  the  normal  equations  for  the  Gauss-Newton  linearization. 
Solving  these  equations  for  g^  minimizes  a  linearized  performance  index 
given  by 

PI2(9i)  =  |lJ(x.)  9i  -  (-f(x.))  |  |2.  (3-11) 

These  normal  equations  (3-10)  are  solved  by  applying  a  Householder 
transformation  to  the  squared  error  shown  in  equation  (3-11).  This 
numerically  more  accurate  approach ®  is  developed  by  introducing  T^T 
into  equation  (3-11)  as  follows: 

PI2(9i)  -  (J  (Xi)  gt  ♦  f(x.))T  TTT( J(x ■ )  9l  ♦  f(xi)) 

PI2(9j )  »  (T  J(Xi)  9i  ♦  T  f ( x  ^ )  )T ( T  J(Xi)  g.  ♦  T  f (Xi ) ) 

PI2(9i)  -  1  IT  J  (xi)  g/Tf  (Xi)  !|2  (3-12 ) 


1 
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where  T  is  the  Householder  transformation.  Because  this  is  an  orthogonal  | 

transformation,  TTT  =  I  and  inserting  this  into  equation  (3-11)  leaves  * 

PI2  unchanged.  In  addition,  T  is  selected  such  that 

T  J  =  (-§-)  (3-13) 

where  R  is  a  5  by  5  upper  triangular  matrix  and  0  is  a  (K-5)  by  5  zero 
matrix.  This  gives 

R  7  2 

PI2  =  IH-§-)  9i  +  (— ~)  II  (3-14) 

where 

T  f  =  (-»)•  (3-15) 

The  gi  which  minimizes  P I 2  is  given  by 

9,  =  -R_1  T  (3-16) 

and 

min  PI2  =  |  |e||2.  (3-17) 

The  details  associated  with  defining  the  elements  of  this  transformation 
T  are  given  in  references  3,  6,  and  7. 

In  the  vicinity  of  a  local  minimum  at  x.+p 

f(xi+1)  =  f(x1  +  g.)  =  f(x.)  +  J(x.)  g.  (3-18) 
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and,  consequently,  (3-11)  becomes 


PI2(gi)  »  Hf(xi+1)U2  =  PIl(xi+1).  (3-19) 

Equation  (3-19)  is  not  a  sufficient  condition  for  a  local  minimum,  but  it 
is  convenient  to  use  it  as  the  termination  criterion  for  the  experiments. 

The  modified  Gauss-Newton  method  is  given  by 

xi+l  =  xi  +  ai  gi  (3-20) 

where  the  positive  scalar  ai  is  selected  such  that 

PIl(xi+1)  <  PIl(x.).  (3-21) 

In  this  method  it  is  said  that  gives  the  direction  of  the  step  while 
ai  gives  its  length. 

Once  g.  is  found,  by  the  Gauss-Newton  method,  PIl(x^)  =  PIl(x. 

+  a.  g.)  is  only  a  function  of  the  scalar  a.;  therefore,  this 
function  can  be  minimized  by  the  selection  of  a^.  There  are  several 
schemes  for  selecting  a-.  The  one  used  in  this  report  finds  the 
minimum  of  a  quadratic  polynomial  which  fits  the  (a^,  PI2)  data  in  the 
g.  direction.  This  is  generated  by  computing  PI2  for  three  different 
values  of  a..  An  analytical  expression  for  the  minimum  in  terms  of 
these  three  data  points  is  derived  in  appendix  8. 

Inspection  of  equation  (3-10)  shows  that  a  non-singular  JTJ  matrix 
is  required  by  the  Gauss-Newton  method.  Ill-conditioning  of  this  matrix 
is  also  a  concern.  The  levenberg-Marquardt  method  is  useful  in  these 
situations,  and  it  also  provides  a  way  to  incorporate  a  priori  informa¬ 
tion  into  the  problem. 


The  Levenberg-Marquardt  method  is  defined  by  equation  (3-5)  but  with 


V  2PIl(xi )  =  2  JT(Xi)  J(x.)  ♦  2  xi  D  (3-22) 

where  a.  is  a  positive  scalar  and  D  is  a  positive  definite  diagonal 
matrix.  By  comparing  equation  (3-22)  with  equation  (3-8)  it  is  seen  that 
2  a.  D  may  be  viewed  as  an  approximation  to  the  term  that  is  omitted  by 
the  Gauss-Newton  method.  The  details  of  the  selection  procedure  for  both 
a.  and  D  are  given  in  reference  4,  but  inspection  of  equation  (3-22) 
shows  that  a ^  =  0  yields  the  Gauss-Newton  method  while  A^  *  ®  with  D 
■  I  yields  steepest  descent. 

In  the  experiments  reported  in  the  next  section,  no  difficulties  with 
the  JTJ  matrix  occurred,  so  the  Levenberg-Marquardt  method  did  not 
prove  useful.  These  experimental  results  are  not  shown  in  this  report. 

Tne  Levenberg-Marquardt  method  is  discussed  because  the  omniu .recti onal 
linear  array  problem  is  singular  prior  to  an  observer  maneuver,  and 
estimates  during  this  time  interval  may  be  required. 
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4.  EXPERIMENTAL  RESULTS 


The  results  from  a  digital  computer  simulation  of  three  different 
sets  of  experiments  are  presented  in  this  section.  The  different  sets 
are  produced  by  changing  the  target's  initial  position.  Each  of  these 
sets  is  composed  of  three  members,  which  are  produced  by  changing  the 
standard  deviation  of  the  measurement  noise.  The  target  parameters  for 
these  three  sets  of  experiments  are  given  in  table  1. 


Table  1.  True  Target  Parameters  with  Speed  in  Meters/Second 
and  Range  in  Meters 


Experiment 

Set 

Vxt 

Vyt 

Rxt(O) 

Ryt(O) 

Rzt(O) 

No.  1 

-5 

10 

10,000 

10,000 

1000 

No.  2 

-5 

10 

40,000 

40,000 

0 

No.  3 

-5 

10 

10,000 

10,000 

10,000 

Experiment  set  No.  1  is  a  short-to-mid  range  situation  while  No.  2  is  a 
mid-to-long  range  one.  Experiment  set  No.  3  may  be  viewed  as  a  single¬ 
bottom-bounce  problem  where  the  actual  target  depth  is  located  by  knowing 
the  bottom  depth.  This  is  shown  in  figure  5. 

The  standard  deviation  for  the  white  Gaussian  measurement  noise  is 
established  from  equation  (2-5)  with  c  =  1500  meters/second  and  b  *  3/8 
meter  (which  is  a  quarter  wave-length  at  a  1-kHz  frequency).  The  three 


standard  deviations  used  in  each  of  the  sets  shown  in  table  1  are 
ot  =.  0.0  seconds,  ot  =  0.05  x  10"^  seconds,  and  ot  =  0.25  x 
10'^  seconds.  The  first  of  these  is  the  noiseless  case  while  the  last 
corresponds  to  a  ^  =  0.01  radians  (approximately  0.5  degrees).  The 
remaining  a t  represents  a  low-noise  ca«e. 

In  all  experiments  the  observer  starts  from  the  origin  at  t  =  0 
seconds  with  a  speed  of  10  meters/ second  and  a  heading  of  0  degrees. 

Just  prior  to  the  measurement  at  240  seconds  the  observer  instantaneously 
changes  heading  to  90  degrees  and  just  prior  to  480  seconds  changes  back 
to  0  degrees.  Each  of  these  headings  defines  a  leg  of  the  problem,  and 
each  leg  consists  of  240  measurements.  The  first  measurement  occurs  at  0 
seconds  while  the  last  occurs  at  719  seconds.  With  the  exception  of  the 
target's  initial  position,  the  target  and  observer  tracks  are  the  same 
for  all  experiments.  Figure  6  shows  these  tracks  projected  onto  the 
(Ry,  Rx)  plane. 

The  time-delay  measurements  are  taken  at  a  1-second  data  rate;  they 
are  then  compressed  by  a  measurement  preprocessor,  which  also  reduces  the 
noise.  For  example,  this  preprocessor  produces  a  "measurement"  at  9.5 
seconds  by  first  adding  together  the  measurements  from  0  to  19  seconds 
and  then  dividing  this  sum  by  20.  These  preprocessed  measurements  are 
then  fed  to  the  estimation  algorithm.  Hereafter,  the  indices  of  all  the 
variables  denote  the  preprocessed  measurement  count;  e.g.,  tm (0 )  is  the 
first  of  these  measurements  and  it  is  associated  with  a  real  time  of  9.5 
seconds  while  tm(l)  is  the  second,  and  real  time  is  29.5  seconds.  The 
error  introduced  by  measurement  preprocessing  is  not  modeled  in  these 
experiments.  This  error  is  zero  only  if  the  true  time  delay  is  constant 
for  each  of  the  20-second  intervals.  Finally,  the  time-delay  pre¬ 
processor  averages  cos  a,  not  a  itself. 
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Figure  5.  A  Single  8ottom-Bounce  Propagation  Path 


Figure  6.  Target  and  Observer  Tracks  in  the  (Ry,  Rx)  Plane 


To  start  the  iterative  least-squares  algorithm  an  initial  estimate  of 
the  target  parameters  is  required.  This  is  specified  by  using  the  first 
measurement  together  with  Vxt  =.  0,  V  t  -  0,  R2t(0)  -  3,000  meters, 
and  Rs(0)  =  20,000  meters;  this  gives, 

Rn(0)  *  c  Rs{0)  tm(0)/b  (4-1) 

and 

Rt(0)  «  *  (Rs2(0)  -  Rzt2(0)  -  Rn2(0))1/2  (4“2) 

where  the  sign  on  R^-(O)  determines  the  initial  side  of  the  cone's 
surface  on  which  the  target  is  assumed  to  by  lying.  Both  signs  are  tried 
and  the  one  which  yields  the  smaller  P 1 1  is  selected  as  the  initial 
point.  Although  this  does  require  that  the  estimate  of  the  target- 
observer  track  be  generated  twice,  it  does  not  require  a  Householder 
transformation. 

Another  initialization  scheme  is  to  select  RT ( 0 )  *  0  and 

Rzt(0)  -  *  (Rs2(0)  -  Rn(0))1/2,  (4-3) 

where  either  sign  may  be  used  without  penalty.  This  scheme  may  also  be 
used  whenever  the  quantity  under  the  square-root  operator  in  equation 
(4-2)  is  negative.  A  future  report  will  investigate  the  use  of  equation 
(4-3)  ana  Ry(0)  *  0  as  the  primary  initialization  procedure. 

In  these  experiments  the  iterative  procedure  is  terminated  when 
PU(x.+1)  is  within  *  10  percent  of  P 12 ( g^ ) .  This  criterion  is 
selected  as  a  convenience  for  these  initial  experiments.  It  is 
recommended  that  a  similar  bound  be  placed  upon  the  component  of  x.+^, 
which  undergoes  the  largest  percent  change.  Monitoring  the  percent 
change  in  PI1  is  also  recommended. 


The  first  experiments  use  the  modified  Gauss-Newton  method.  The 
initialization  uses  the  sign  on  Ry(O)  which  produces  the  smaller  PI1, 
so  the  estimates  associated  with  generating  the  additional  target- 
observer  track  are  not  shown  in  the  results.  However,  some  additional 
experiments  do  show  the  effects  of  using  the  "wrong"  sign  on  Ry(0). 

This  is  followed  by  some  results  for  both  a  one-legged  and  a  two-legged 
problem.  Finally,  results  obtained  with  the  Gauss-Newton  method  are 
mentioned . 

The  results  from  the  first  experiment  set  are  shown  in  tables  2a,  2b, 
and  2c;  from  the  second  in  tables  3a,  3b,  and  3c;  and  from  the  third  in 
tables  4a,  4b,  and  4c.  In  all  of  these  tables  the  iteration  number 
defines  the  number  of  calls  to  the  Householder  transformation  subroutine, 
and  the  last  row  shows  the  PI1  which  occurs  when  the  true  target  para¬ 
meters  are  used.  As  expected,  this  PI1  is  larger  than  the  one  reached  by 
the  algorithm  in  all  cases. 

Inspection  of  all  these  tables  shows  that  the  algorithm  converges  by 
the  fourth  iteration,  except  for  the  noiseless  cases.  In  these  cases  it 
is  seen  that  Rzt(0)  is  the  only  target  parameter  wh ich  experiences  more 
than  a  10  percent  decrease  in  moving  from  the  third  to  the  fourth 
iteration,  and  this  only  occurs  in  tables  2a  and  3a.  Perhaps,  after  the 
fourth  iteration  a  1-dimensional  search  or  a  less  stringent  performance 
bound  should  be  used  for  R  ^(0).  Regardless,  further  investigation  or 
discussion  of  these  noiseless  cases  is  not  warranted. 

Inspection  of  tables  2a,  2b,  and  2c  shows  that  the  magnitude  of  the 
error  in  the  final  estimate  increases,  as  expected,  when  the  variance  of 
the  noise  is  increased.  This  is  also  shown  in  the  other  sets  of  tables. 

A  Monte  Carlo  simulation  or  an  analytical  study  could  quantify  this 
observation . 
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Table  2a.  Results  from  the  Modified  Gauss-Newton  Method 
with  ct  «  0.0  Seconds  (First  Member  of  Set  No.  1) 


Iteration 

Number 

Vxt 

vyt 

Rxt(0) 

Ryt(0) 

R2t(0) 

PI1 

PI2 

0 

0 

0 

13,800 

14,161 

3000 

3.79X10® 

4.70X105 

1 

2.75 

10.25 

6179 

6993 

3087 

1.89X107 

8463 

2 

-2.80 

12.99 

8742 

8563 

1701 

5.87X10® 

46.47 

3 

-5.39 

9.64 

10,230 

10,233 

1127 

3561 

0.024 

4 

-4.99 

10.01 

9997 

9997 

1006 

4.58 

0.057 

5 

-5.00 

10.00 

10,002 

10,002 

1002 

0.059 

— 

True 

Parameters 

-5.00 

10.00 

10,000 

10,000 

1000 

0.57 

_ 

Table  2b.  Results  from  the  Modified  Gauss-Newton  Method  with 
ot  =  0.05  X  10_5  Seconds  (Second  Member  of  Set  No.  1  ) 


Iteration 

Number 

Vxt 

vyt 

Rxt(0) 

Ryt(0) 

Rzt(0) 

PU 

PI2 

0 

0 

0 

13,827 

14,136 

3000 

3.83X10® 

4.75X105 

1 

2.70 

10.25 

6178 

6993 

3109 

1.86X107 

8096 

2 

-3.05 

12.89 

8899 

8699 

1734 

5.88X105 

763 

3 

-5.56 

9.53 

10,328 

10,326 

1165 

4014 

734 

4 

-5.18 

9.89 

10,106 

10,100 

1054 

737 

— 

True 

Parameters 

-5.00 

10.00 

10,000 

10,000 

1000 

1029 

____ 

Table  2c.  Results  from  the  Modified  Gauss-Newton  Method  with 
ot  ■  0.25  X  10"5  Seconds  (Third  Member  of  Set  No.  1) 


Iteration 

Number 


Vxt 


0 


Rxt<0) 

Ryt(O) 

Rzt(O) 

PI1 

13,801 

14,161 

3000 

3.79X108 

6285 

7058 

2939 

2.38X107 

7637 

7572 

1526 

6.52X105 

9494 

9451 

858 

37,858 

9191 

9166 

613 

29,739 

10,000 

10,000 

1000 

37,793 

PI2 


5. 02X10& 


42,580 


29,736 


29,648 


Table  3a.  Results  from  the  Modifiea  Gauss-Newton  Method  with 
at  *  0.0  Seconds  (First  Member  of  Set  No.  2) 


Iteration 

Number 

Vxt 

0 

0 

1 

17.56 

2 

-5.74 

3 

-4.85 

4 

-5.02 

5 

-5.00 

6 

-5.00 

7 

-5.00 

8 

-5.00 

True 

Parameters 

-5.00 

Rxt(O) 


13,812 


Ryt(0 


14,150 


26,916 


41,203 


39,766 


40,032 


39,996 


40,003 


40,004 


40,004 


40,000 


Rzt(0) 


3000 


4659 


2628 


1013 


454 


178 


101 


89 


87 


Pli 


5.51X107 


3.64X106 


1.37X105 


3813 


98 


3 


0  .173 


0.163 


0.134 


missing 


PI  2 


22,182 


2041 


0.148 


0.148 


0.140 


0.142 


0.137 


0.134 


Table  3b.  Results  from  the  Modified  Gauss-Newton  Method  with 
at  =  0.05  X  10“5  Seconds  (Second  Member  of  Set  No.  2) 


Iteration 

Number 

Vxt 

V 

Rxt(O) 

Ryt(O) 

Rzt(0) 

PI1 

P 12 

0 

0 

0 

13,806 

14,156 

3000 

5.46X10? 

30,848 

1 

17.60 

30.94 

26,416 

26,921 

4668 

3.62X106 

25,604 

2 

-5.70 

9.92 

41,177 

41,016 

2663 

1.54X105 

23,343 

3 

-4.83 

9.92 

39,558 

39,616 

1131 

26,524 

23,343 

4 

-4.98 

9.91 

39,820 

39,844 

754 

23,378 

— 

True 

Parameters 

-5.00 

10.00 

40,000 

40,000 

0.0 

24,205 

— 

Table  3c.  Results  from  the  Modified  Gauss-Newton  Method  with 
ot  -  0.25  X  10-5  Seconds  (Third  Member  of  Set  No.  2) 


Iteration 

Number 

Vxt 

vyt 

Rxt(O) 

Ryt(O) 

Rzt(O) 

PI1 

PI2 

0 

0 

0 

13,757 

14,204 

3000 

5.04X107 

2.03X105 

1 

19.54 

32.84 

25,919 

26,460 

4549 

3.83X106 

3.31X105 

2 

-7.51 

7.22 

41,133 

41,199 

2844 

4.00X106 

3.04X105 

3 

-6.94 

7.14 

40,120 

40,346 

1880 

3.06X106 

— 

True 

Parameters 

_ 

-5.00 

10.00 

40,000 

40,000 

0.0 

3.77X105 

_ . _ 

) 
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Table  4a.  Results  from  the  Modified  Gauss-Newton  Method  with 
ot  *  0.0  Seconds  (First  Member  of  Set  No.  3) 


Iteration 

Number 

Vxt 

v>t 

Rxt(O) 

Ryt(O) 

Rzt(O) 

PI1 

PI2 

0 

0 

0 

16,039 

11,565 

3000 

8.69X108 

3.91X10& 

1 

-1.66 

1.58 

13,726 

11,045 

10,376 

4.63X108 

1.49X108 

2 

-6.58 

6.73 

9975 

10,586 

10,905 

1.34X107 

1163 

3 

-5.16 

10.11 

10,196 

10,112 

10,059 

47,132 

0.0066 

4 

-5.00 

10.00 

10,000 

10,000 

10,001 

0.522 

0.0073 

5 

-5.00 

10.00 

9999 

10,000 

10,000 

0.0061 

0.0061 

6 

-5.00 

10.00 

9999 

10,000 

10,000 

0.0060 

— 

True 

Parameters 

-5.00 

10.00 

10,000 

10,000 

10,000 

0.261 

Table  4b.  Results  from  the  Modified  Gauss-Newton  Method  with 
at  »  0.05  X  10~5  Seconas  (Secono  Member  of  Set  No.  3) 


Iteration 

Number 

Vxt 

vyt 

Rxt(O) 

Ryt(O) 

Rzt(O) 

PI1 

P 12 

0 

0 

0 

16,036 

11,569 

3000 

8.68X108 

3.80X108 

1 

-1.67 

1.58 

13,730 

11,049 

10,375 

4.63X10® 

1.42X105 

2 

-6.55 

6.75 

9950 

10,558 

10,886 

1 . 39X1 07 

2579 

3 

-5.00 

10.12 

10,083 

10,017 

9980 

49,698 

1815 

4 

-4.86 

10.01 

9900 

9915 

9931 

1816 

— 

True 

Parameters 

-5.00 

10.00 

10,000 

10,000 

10,000 

2003 

— 
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Table  4c.  Results  from  the  Modified  Gauss-Newton  Method  with 
<j£  «  0.25  X  10~5  Seconds  (Third  Member  of  Set  No.  3) 


Iteration 

Number 

m 

Vyt 

Rxt(O) 

Ryt(O) 

Rzt(0) 

PI1 

PI  2 

0 

0 

0 

16,065 

11,529 

3000 

8.77x10s 

4.92X105 

I 

-1.70 

1.54 

13,718 

10,995 

10,364 

4.65X10® 

2.28X105 

2 

-6.87 

6.54 

10,059 

10,609 

10,890 

1.44X107 

55,891 

3 

-5.67 

9.97 

10,440 

10,247 

10,152 

1.01X105 

53,939 

4 

-5.46 

9.86 

10,213 

10,112 

10,072 

53,937 

— 

True 

Parameters 

-5.00 

10.00 

10,000 

10,000 

10,000 

60,327 

— 

Comparing  tables  2c,  3c,  and  4c  shows  that  the  final  value  of  PI1 
increases  with  increasing  slant  range.  The  results  shown  in  these  tables 
demonstrate  that  the  modified  Gauss-Newton  method  converges  to  a  reason¬ 
able  estimate  of  the  target  parameters  in  approximately  four  iterations. 

The  high  noise  cases  for  these  three  experiment  sets  were  repeated, 
but  they  were  initialized  on  the  "wrong"  side  of  the  cone;  i.e.,  the  side 
which  yields  the  larger  PH.  The  one  from  experiment  set  no.  3  displayed 
results  similar  to  those  reported  in  table  4c  while  the  remaining  two 
high  noise  cases  performed  poorly.  After  five  iterations  the  one  from 
experiment  set  no.  2  was  stopped  because  convergence  to  a  reasonable 
estimate  was  judged  unlikely  After  seven  iterations  the  one  from 
experiment  set  no.  1  was  also  stopped;  however,  convergence  in  one  more 
iteration  was  expected  but,  unfortunately,  with  the  wrong  sign  on  the 
Rzt(0)  estimate. 
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To  investigate  the  effects  of  fewer  legs,  the  third  experiment  set 
shown  in  table  1  was  repeated,  first  with  two  legs  (0  to  479  seconds)  and 
then  with  one  (0  to  239  seconds).  Table  5  shows  the  final  results  of  the 
two-legged  problem. 


Table  5.  Estimates  from  the  First  Two  Legs  of  Experiment  Set  No.  3 


Number  of 
Iterations 

at 

Vxt 

vyt 

Rxt(O) 

Ryt(O) 

RZt(0) 

PI1 

5 

0.0 

-5.0 

10.0 

10,000 

10,000 

10,000 

0.0033 

4 

0.05X10-5 

-6.2 

9.50 

10,820 

10,726 

10,629 

965.0 

7 

0.25X10-5 

10.0 

2.55 

-2495 

_ 

1773 

-90 

852.0 

Tnis  table  shows  that  the  noiseless  and  the  low  noise  cases  converged  to 
reasonable  estimates  while  the  high  noise  case  converged  to  a  worthless 
result.  This  worthless  result  may  be  caused  by  fewer  legs  or  by  fewer 
measurements;  additional  experiments  could  resolve  this  point.  As 
expected,  the  one-legged  problem  generated  diverging  results  because  the 

matrix  is  singular.  The  Levenberg-Marquardt  method  could  be  used 
if  an  estimate  on  this  leg  is  required. 

Finally,  the  Gauss-Newton  method  was  used  to  treat  the  three  high 
noise  cases.  The  results  obtained  were  essentially  the  same  as  those 
shown  in  tables  2c,  3c,  and  4c;  consequently,  the  additional  computations 
required  to  select  a  step  size  are  not  justified  by  these  experiments. 
However,  these  experiments  do  not  constitute  an  exhaustive  set  of  all 
possible  target-observer  situations,  and  the  modified  Gauss-Newton  method 
may  yet  prove  useful. 
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5.  CONCLUSIONS  AND  FURTHER  WORK 

The  main  conclusions  are: 

1.  There  is  a  sign  ambiguity  on  the  estimate  of  the  target's 
depth . 

2.  The  problem  is  singular  on  the  first  leg. 

3.  Both  the  Gauss-Newton  method  and  the  modified  Gauss-Newton 
method  converge  in  four  iterations  for  these  experiments. 

4.  The  Jacobian  matrix  is  singular  if  the  estimate  of  the 
target's  depth  is  zero. 

5.  Initialization  of  the  iterative  estimator  affects 
convergence. 

The  sign  on  the  target's  depth,  *^zt(0) ,  cannot  be  resolved  by  an 
observer  operating  at  a  constant  depth.  This  is  seen  from  equation 
(2-16),  and  it  was  also  demonstrated  experimentally  although  these 
results  are  not  shown. 

The  target  tracking  problem  with  a  linear  array  of  omnidirectional 
elements  is  singular  on  the  first  leg.  This  was  shown  in  the  experi¬ 
ments.  An  analytic  treatment  of  this  point  is  not  a  high  priority  item, 
but  is  recommended  since  it  would  lend  additional  insight  into  this 
tracking  problem. 

Both  the  Gauss-Newton  and  the  modified  Gauss-Newton  methods  demon¬ 
strate  convergence  in  approximately  four  iterations.  Additional 
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experiments  are  needed  to  investigate  convergence  and,  in  those 
situations  where  convergence  occurs,  to  determine  if  convergence  in  four 
iterations  is  typical. 

Inspection  of  appendix  A  shows  that  the  last  column  of  the  Jacobian 
matrix  is  zero  if  the  estimate  of  the  target's  depth  is  zero;  i.e.,  the 
target  is  estimated  to  be  at  the  same  depth  as  the  observer.  This  sing¬ 
ularity  will  be  investigated  in  a  subsequent  report. 

Initialization  of  the  estimation  algorithm  can  reduce  the  number  of 
iterations  needed  for  convergence.  Further  experiments  on  this  item  are 
needed. 

To  evaluate  the  Gauss-Newton  method,  more  experiments  with  different 
target-observer  tracks  and  with  different  noise  variances  are  needed.  A 
Monte  Carlo  simulation  is  also  needed. 

Oifferent  PI1  functions  should  be  tried,  with  a  view  towards  finding 
one  that  reduces  the  bias  in  the  final  estimate.  Also,  a  probabilistic 
description  of  this  problem  might  be  more  useful  than  the  least-squares 
approach  used  here. 

Experiments  which  expand  the  iterative  least-squares  estimator  to 
inc.uded  measurements  from  additional  sensors  would  be  useful. 

Finally,  results  from  an  extended  Kalman  filter  which  explores  the 
effects  of  different  signal  models  should  be  documented. 
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APPENDIX  A.  THE  JACOBIAN  OF  f(x) 

The  Jacobian  of  f(x),  which  is  used  in  the  performance  index  of 
equation  (3-1),  is  shown  below.  This  f(x)  is  defined  by  equation  (3-2). 

The  Jacobian  is  defined  by 
J  =  (3Vaxj) 

for  j  =  1,2, 3, 4,  and  5;  and  k  =  0,1 ,2, . . .  ,35. 

The  entries  in  the  kth  row  of  J  are  given  by: 

3fk/ax1  =  3VaVxt  *  kh  3  VsRxt(°) 

3fk/3x2  .  3VaVyt  =  kh  3  fk/3Ryt*0^ 

3fk/sx3  =  afk/ aRxt (0)  * -sin  (k)  -  ctm(k)  Rx(k)/b  R$(k) 

3fk/3x4  =  3fk/3Ryt(0)  = -cos  *  (k)  -  ctm(k)  Ry(k)/b  R$(k) 

3fk/3x5  =  3fr/3Rzt(0)  =  -ctm(k)  Rz(0)/b  R$(k). 
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APPENDIX  B.  STEP  SIZE  FOR  THE  MODIFIED  GAUSS-NEWTON  METHOD 


In  this  appendix  a  scheme  is  given  for  selecting  the  step  si2e  a^ 
in  the  modified  Gauss-Newton  iterative  formula. 


i+1 


-  xi  +  ai 


9i 


(B-l ) 


Actually,  because  it  is  not  normalized,  the  direction  g^  also  con¬ 
tributes  to  the  size  of  the  step.  It  is  convenient  to  redefine  equation 
(B-l)  as 


where  a,  denotes  the  jth  value  of  the  step  size  at  the  ith  iteration. 

J 

Once  gi  is  found  from  the  Gauss-Newton  equations,  PI1  is  a  function 
only  of  a^;  i .e. , 

PIKaj)  «  PIl(x1  +  9j  gi)  ( B— 3 ) 


and  this  is  minimized  by  a  judicious  selection  of  a..  In  this  report, 

J 

a^  is  defined  by  the  minimum  of  a  quadratic  polynomial  which  passes 
through  three  (a,,  PIl(a.)}  data  points.  For  equally  spaced  values 

J  J 

of  a.,  the  minimum  of  this  quadratic  is  given  by 
J 

(a2+a3)  Pll(a^)  -  2(a1+a3)PIl(a2)  +  (a^)  PI1  (a3 ) 
am  "  2  P 1 1  ( a1 )  -  4  Pll(a2)  +  2  PTITSJJ 

where  a3  >  a2  >  a^  1  0* 


The  first  of  these  data  points  is  readily  available,  namely,  (a^  * 
0,  PI1(0)  *  P 1 1 ( x i ) );  ano  if 


PI1(1)  <  PI1(0), 


( B— 5 ) 


a„  -  1  gives  the  second  data  point  and  a3  =  2ag  gives  the  third. 

However,  if  equation  (B-5)  is  not  satisfied  the  length  of  the  interval  is 

reduced  by  selecting  a2  =  1/2  and  a3  =  2a2  =  1,  provided 

PIl(l/2)  <  PI1(0).  (B'6 

If  this  is  unsuccessful,  the  next  selection  is  a2  =  1/4  and  a3  = 

2a?  =  1/2,  and  subsequent  selections  are  given  by  repeatedly  dividing 
a2  by  2.  This  continues  until  P I 1 ( a£ )  <  Pll(a^)  or  a  threshold 
which  causes  termination  of  the  estimation  algorithm  is  crossed.  After 

is  found,  then  PIl(am),  PIl(a2) ,  and  PIl(a3)  are  compared  to 
determine  which  of  these  is  the  smallest.  This  is  necessary  because  the 
quadratic  polynomial  may  not  always  provide  a  good  fit  to  the  PI1 
function  and  PIl(a2)  or  PIl(a3)  may  be  smaller  than  PIl(am).  A 
diagram  of  this  scheme  is  shown  in  figure  B-l. 


APPENDIX  C.  NUM8ER  OF  OPERATIONS  FOR  THE  MODIFIED  GAUSS-NEWTON  METHOD 

The  number  of  operations  (additions,  multiplications,  divisions, 
square  roots,  and  trigonometric  functions)  required  for  one  iteration  of 
the  modified  Gauss-Newton  method  is  determined  in  this  appendix. 

This  method  is  defined  by 

xi+l  *  x.  ♦  a,  g,  (C-l) 

where  a^  is  the  step  size  and  is  determined  from 

J(*1)  91  -  f(Xi)  (C— 2 ) 

by  applying  a  Householder  transformation  to  the  (J(x.)  j  f(x-)) 
array.  This  transformation  gives 


R  9i  -  f  (C-3) 

where  R  is  an  upper  triangular  matrix  and  g^  is  found  by  back 
substitution. 

The  number  of  operations  needed  to  determine  ai  is  given, 
primarily,  by  the  number  of  additional  times  that  PI1  must  be  evaluated. 
Appendix  B  showed  that  a  minimum  of  two  additional  evaluations  per 
iteration  are  required,  and  the  experiments  demonstrated  that  frequently 
this  was  all  that  was  needed.  Furthermore,  equation  ( B— 4 )  shows  that  in 
addition  to  these  operations,  5  additions,  6  multiplications,  and  1 
division  are  required.  The  divisions  needed  to  reduce  the  step  size  (see 
figure  B— 1 ) ,  assuming  that  a  reduction  is  called  for,  are  not  counted. 
Therefore,  the  minimum  number  of  operations  required  to  find  the  step 
size  are  2(K-l)+5  additions,  2  K+6  multiplications,  and  1  division. 
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The  equations  which  define  the  elements  in  the  kth  row  of  the 
(J(x.)  j  f(x.))  array  are  given  in  appendix  A.  The  operations  needed 
to  evaluate  these  elements  are  shown  in  table  C-l.  The  operations  needed 
for  PI1  are  included. 


Table  C-l.  Number  of  Operations  Required  for  the  ( J  If )  Array  and  P 1 1 . 


Equation 

4- 

X 

# 

r~ 

trig 

f (k)  =  Ru(k)  -  R$(k)tj(k)c/b 

1 

2 

- 

- 

- 

R11(k)  *  -Rx(k)sin*(k)  -  R^(k)cosiJi(k) 

1 

2 

- 

- 

2* 

Rs(k)  -  (Rx2(k)  ♦  Ry2(k)  ♦  Rz2(0))1/2 

2 

3 

- 

i 

- 

af/3Rxt(°)  =-sim|i(k)  -  RJ((k)t1  (k)c/R$(k)b 

1 

2 

1 

- 

*★ 

af/sR^t(0)  =-cos*(k)  -  R  (k)t1(k)c/Rs(k)b 

1 

1 

★★ 

- 

*★ 

3f/aRzt(0)  =  -  Rz(0)t1(k)c/Rs(k)b 

- 

1 

★  * 

- 

- 

af/aVxt  =  kh  af/a  Rxt(0) 

- 

2 

- 

- 

- 

af/avyt  =  kh  af/a  Ryt(0) 

- 

1 

- 

- 

- 

Total  for  one  row 

6 

14 

1 

i 

Total  for  K  rows 

6K 

14K 

K 

K 

★ 

P 1 1  requires 

K- 

1  K 

- 

- 

- 

♦Only  2  per  leg,  so  they  are  not  included 

♦♦c/b  is  computed  off-line,  so  c/rs(k)b  is 
tl (k)c/R${k)b  is  computed  only  once. 

in 

1 

the  totals. 

division.  Also 

The  number  of  operations  needed  to  find  R  (k)  and  R  ( k ) ,  given  an 

x  y 

estimate  of  the  target's  velocity  and  initial  position  and  the  observer 
track,  are  not  counted  because  they  are  the  same  for  all  iterative 
estimation  algorithms  which  use  a  Cartesian  signal  model. 
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The  Householder  transformation  used  in  the  iterative  least-squares 
algorithm  is,  in  fact,  a  product  of  Householder  transformations.  The 
first  member  of  this  product  operates  on  the  matrix 

A  =  (J  j  f),  (C— 4 ) 

which  has  K  rows  and  n  columns,  where  n  ■  6.  The  second  member  in  this 
product  leaves  the  first  row  and  column  of  the  transformed  A  matrix 
unchanged,  so  it  operates  on  a  (K— 1 )  by  ( n— 1 )  matrix.  The  next  member 
operates  on  a  (K-2)  by  (n-2)  matrix  and  so  on  until  the  first  ( n— 1 ) 
columns  of  A  are  transformed.  The  number  of  operations  needed  to 
transform  A  so  that  its  new  first  column  has  the  required  structure  is 
shown  in  the  following  table.  The  symbols  used  in  this  table  are 
consistent  with  the  Householder  transformation  algorithm  shown  on  page 
63  of  reference  3,  except  that  K  is  used  here  instead  of  m. 


Table  C-2.  Number  of  Operations  Required  to  Transform  the  First 
Column  of  the  A  Matrix 


Symbol 

+ 

X 

l 

r~ 

s 

K-l 

K-l 

- 

1 

u(l) 

1 

- 

- 

- 

B 

- 

1 

1 

- 

y 

(K-l)  (n-1) 

(K*l)  (n-1) 

- 

- 

A 

K(n-l) 

K(n-l) 

- 

- 

Total 

Kn+(K-1 ) (n-1 ) 

2Kn-*-n-K+l 

1 

1 

In 

the  omnidirectional  least- 

■squares  problem, 

n  =  6,  hence. 

the  first 

five  columns  of  A  must  be  transformed  to  upper  triangular  form.  The  next 
table  shows  the  number  of  operations  needed  for  each  of  these 
transformations. 


i 


< 

i 


Table  C-3.  Number  of  Operations  Required  for  the 
Householder  Transformation  of  A 


Transformation 
Defined  by 

Column  No. 

+ 

X 

r~ 

1 

Kn+(K-l)(n-l) 

2Kn+n-K+l 

1 

i 

2 

(K-l)  (n-1 )+(K-2 ) 

(n-2) 

2 (K-l)  (n-1 )+n-K+l 

1 

i 

3 

(K-2)  {n-2 )+(K-3 ) 

(n-3) 

2 (K-2)  (n-2)+n-K+l 

1 

i 

4 

(K-3)  (n-3 )+(K-4 j 

( n— 4 ) 

2 (K-3)  (n-3)+n-K+l 

1 

i 

5 

(K-4)  (n-4)+(K-5) 

(n-5) 

2 ( K-4 )  (n-4 )+n-K+l 

1 

i 

Tota  1 

10Kn-25K-25n+85 

10Kn-25K-15n+65 

5 

5 

Total  with  n  »  6 

35K-65 

35K-25 

5 

5 

To  compute  PI2  requires  ( K— 6 )  additions  and  (K-5)  multiplications. 

The  back  substitution  algorithm  requires  10  addtions,  10 
multiplications,  and  5  divisions. 

The  total  number  of  operations  required  for  one  iteration  of  the 
modified  Gauss-Newton  method  is  shown  in  table  C-4. 


Table  C-4. 

Final  Count  of  Operations  for  One 
Modified  Gauss-Newton  Method 

Iterations 

of  the 

Function 

♦ 

X 

i 

r~ 

Step  size 

2K*3* 

2K+6* 

1* 

- 

(J  j  f) 

6K 

14K 

K 

K 

PI1 

K-l 

35K-25 

_ 

* 

♦minimum 


j 

1 


C-4 


Table  C-4.  Final  Count  of  Operations  for  One  Iterations  of  the 
Modified  Gauss-Newton  Method  (Cont'd) 


Function 

+ 

X 

* 

J— 

Householder 

Transformation 

35K-65 

35K-25 

5 

5 

PI2 

K-6 

K-5 

- 

- 

Back  substitution 

10 

10 

5 

- 

Gauss-Newton  Total 

43K-62 

51 K— 20 

K+10 

K+5 

Minimum  Number  of 
Extra  Operations 
for  Moaified  G-N 

2K+3 

2K+6 

1 

0 

When  faced  with 

an  ill-conditioned  problem. 

it  is  advisable  to 

inter- 

change  the  columns  (and  maybe  the  rows  as  well)  of  J  before  applying  the 
Householder  transformation.  This  is  explained  in  reference  6.  The  number 
of  operations  required  for  this  interchange  are  now  shown  in  table  C-4. 


C-5/C-6 
Reverse  Blank 


APPENDIX  D.  COMPUTER  PROGRAM  LISTING 


This  appendix  contains  the  computer  programs  used  to  generate  the 
results  shown  in  section  4.  There  are  two  parts  to  this  program:  the 
first  generates  the  averaged  measurements;  the  second  estimates  the  target 
parameters.  Six  significant  digits  are  carried  in  all  computations. 


0  LPRlNT-TnI3  PREGRAM  generates  Trl£  averaged  measurements- 
l  DATA  10.10.iT»U>90.0.-10.U.ia»-5»10»iOOOC»100DC»10CO»C»0 


a  1  TRANSLATION  C7  THE  NAMES  Q F  THE  VA.IAcLES  USED  IN  THIS  PROGRAM 

10  SO-SFEED  OF  THE  OBSERVER .  INPUTTED  AS  SL.SUSZ. 
is  GAMAC-HEADINC  Of  THE  OBSERVER .  IN.-oTTED  AS  CD.G1.G2 

12  VZO=Z  VELOCITY  OF  THE  OS  SERVER.  INF'uTTED  AS  20Q.Z10.Z20 

13  ' VXT-X  vELCCI'Y  OF  THE  TARGET 

n  • uyt ~y  *  ■ 

IS  1 RXT -X  POSITION  OF  THE  TAkGET 

14  1 RYT  *Y 
17  ‘RZT-Z 

ie  'Vixa*x  velocity  of  the  observer 

IV  'VZYO*Y 

::  'Rixo*x  position  of  the  observe* 

21  R2YO-Y 

22  'R3Z0-Z 

23  R'AX-nELATIVI  X  POSITION  (  TARGE'-OE  SERVER) 

2i  'P5Y=  *  Y  * 

20  ‘ R6Z=  •  2 

24  MEAN-MEAN  OF  MEASUREMENT  NOISE 

27  'SOEV-STANCAAD  OEVIATICN  OE  THE  MEASUREMENT  NOISE 
23  1  Cl 1-GAM AC- 161  IN  RADIANS.  ANGLE  OF  ROTATION 

29  'H-SAMFvING  PERIOD 

3)  • RA-R1 1 . »AR AULEL  TO  2ENTEF  Life  OF  ARRAY 

21  '  Rp -R T  »  F  E  VPLpO-CU- AH  - 

22  AlFHA-OOS1 "P-E  ANGLE. 

3 j  1  TOE- AY »TR. Jt  TIME  Ol-mF 

3-  SI* JSCS  IN  RANDOM  N  Jit  Et\  GENERATOR 

so  meise-  •  • 

30  ' TMEAS'NUISY  TIME  DELAY 

37  Av TM-CNE  OCMLnSIGnal  VA.-IA.-lE  .26  ELEMENTS. AVERAGED  TMEAG 

:e  ■  stm-as-used  to  gei  e..a'o  -  avtm 

IOC  *«*  INF  jT  OBSERVE*  Ij.ITIal  .-A.vAMETERS  AND  MANEUVERS.  In-v" 

TAR  52  T  F  ARAMS'  EhS.  Inf  j'  MEAN  AN  D  S  LEV  CF  MEAijHEMEfT  NOISE.*** 

1C3 
i  :m  ■ 

1.5  PRI:.I-DA*A  E  AYE  1 . 2 .  C'C  .  —  32  GE  *vER  Sf  tED.  C  *  l  >  L  —  EES  .  HEADING  .0.1.2 - Z  VE.L 

C'7f  : ,  :,Z _ TA\SE~  >  iy  e  —  , .  f  ■»  El  .  .  f  GS , *  r  FOS..Z  PCS.. - — — NjlSE  mEAn.SClV- 

TY.L  CChI  ...  PER  data  £Nr£;Eu* 


T*TSP*G£  ^  Pr 
fhvw  ’  1 


-7  y:Mvn  mcxicxBi* 

_  *  J  Ju  •* 


D-l 


1C  I: 


**»  ST  Or  STmT'EM  A  _LC.-*S  uS£R  TO  EVER  DATA.  T^ffi  OGn  ' 
TC  C'Cn’InE  OvCE  DATA  IS  EMC  IE*  *«* 

STOP 


1 17 
HE 


125 

133 

121 

132 

133 

1 34 

135 

136 
107 
138 
HC 

142 

143 
145 
143 

149 

150 

151 

152 

153 
155 
15b 
159 

1©C 
161 
162 
163 
195 
1 9  © 
197 
199 
2)0 
238 


'**«  READS  INr-T  DA“h  CC-IN23  If1  REr  lOC.cvc 
**»  INITIAL  ObSEF  vLF  Sf'EEC  A.sC  HEADING  (IN  DEGREES 
ArVC  RADIANS!.  *** 

R 1 4L  3  0  *  S  1  *  S2  f  GO  tGl  »G2»2CC  »2lC»22C»v<T  *  lv  f  T  *  R  XT  »RYT » R.Z7  » F  c.  An  »  £  *.  i 
SO«SO :&A*AC*3C  :n-i»G3*.Ql7A533 


•  *»*  sampling  fef.eid  is  h  seconds.*** 

H-l 

•  ran:om  is  US CO  tY  the  noise  generator 
ram)  Of 

•  ***  INITIAL  X#Y>AN0  Z  OJSERvER  SPEEDS  ARE  SET.  *** 
UlXQ*Sa«SIMNrt>  :92Y0*SU*CSS<wH>  :VZQ=ZCC 

i*t  ELIMINATES  IN IT IA .  QO "SET  IN  cINE  320. *«* 

R 1  :<C=--v i XO «?1 : RZ Y 3  * -v2vJ*h : R 3ZC« -L ZO  «  H 

•  INPUT  FCR  A ‘.ERASED  T*EA3  L50:  W-.Hr.  STARTS  A"  7QC 
OIK  AUTK<3H  J SThEAr«0  Su®-3  SK*0 

THERE  ARE  3©  ELEftEvTS  IN  A‘)TM «  0-EwE-4EnT  C3R  ?  ESPOnDI ''G  41  "h 
7IM£  =  ?  .  5  SCOINl'S  a*  l-  3 1-  Id  ■  Ew  r  Wl~,i  ’C».5  SCI. 

'  ***  ELIMINATES  INI'I A _  QFTIE~  IN  _I'«E  3«C  .  *** 

RXT*-VX“*H*RXT':hY7*-LV?.«r'«.Rr  : 


LPRIN?  'K'}'  •?#R4X'J< 

•  5  *  TOELAY • » ■ TMEAS  * 


* ;  * R5Y  * ;  • 


*  F  6l ' 


•FA' i • 


’»»*  START  OF  MAIN  LOOP.  RUNS  T>R  El  l£G3»  4  HlMdlS  EAIH.*«* 
■  STEP  ON  COUNTER  »;  IS  h.  TmE  3A«F-.INC  PERIOD 


FOR  K«U  TC  723  STEP  1 

’  ***  CHECKS  "CF  FIRST  OEiE'vEK  HANj.vt  •  Am  SE“ 

NEH  SPEED  AND  HEADING. ««* 

237  * 

24  2  I r  24{  ■  *k  ANl  4  9)  7F_N  245  £i  ££  2 cl 
245  SO^-Si  I  )-frtC*Gl  JVZC-ZIC 

IDE  *  *»*  CH  1103  FOR  SECIsC  J:S..A*cr  PnC  S . .  ~  9 

NEW  SPEED  AN)  >’lA  Mn;.  •»* 

16  2  I,"  r.  =A90  TmEn  SI*Si:  GAriA  1  IC*J  Cl 

Z7d  '  »**  rlsCG  ROTATION  AN)_t  "  CR  "hC  ~  .  F I  • -A  .  I  )  *  •  *  * 

«**  r  IMS  DDDE  wc  •:  mIhJInO  IN  h  FO  I  V  ’ 

2  3l  Cl  1  =  <  GArA3-  lc  3  M  .  0  :  J  J  :  DA  ">AI  •. 


0-2 


w5f*SIisK3T^un^u-' 


} 

I 


:'E  '**»  X . v  .  ANC  2  UE.-.  AND  F:zz~z:t  O'  'HE  ;_s  E;. £:.•** 

IXtaxXXXV'CXXl  k  (Ikllrf  Iilxllkkiuklklkklllklltak 


30C  91'<E-SE»SE'-k ,  10  2  Y0=$0xC3£'  GAiA J ' 

32  C  RlXO='mQ*-t.hlxC  :R2YO92Y0*HhR  2v2:f  322*vZ3*-"»F3ZG 
3?  J  1 

333  «**  X.Y.AnO  l  TARGET  FOiITI3N.*«» 

ixxitK(vxx(r<kimixiixxiixxai 

IXmxfliKIIttXISXXIKXXIflkfkXI 

34:  rxt=v<t* i*ax7:.gyt*uyt*.9*ryt:r;t*rz" 

323  '*X»  RELA7I oZ  X.YfANC  2  POSITION  <  TA-.G"-03 SEN /IF ;  .  **» 

llkXXIfk  IXIXIVXXXiXXKXXY 


360  R  AX  sxaT-x  1 XO  *  ft 3  r  »R  YT-R2  YC  2  R  oZ*RZT-R3ZC 
3a  1  ' 

36Z  '  ***  RELATIVE  XiY.Z  Aha  FORnEO.  *** 

*  «»  TFANSFGE  1A7ILN  TO  '0»ED  *** 

■  II  AR-AY  CCCRJIxA'Ea  IS  NEXT.  xr»* 

363 

3a6 

367  • - 

SO:  '•  TRAS.-CR.iAIICN 

523  RA=RA<»SIr»(CU;-R5Y»C3i(Cli) !  R>=R“  <«CC3:C11  >-R5Y*3IM  Oil  > 

5*  )  RSMR.A:tRA*RF«<r  +RoZ*f,6Z )  C  .1 


SI2  '  *«*  TRUE  TIME  CELAY  *** 

260  A-PHAaRA/RS 

sec  T3£LAY=ALF,-.A/100C  . 

332  • - 

63u  ' 

620  '***  GAUSSIAN  NOISE  GEnERA'OR.  *** 

63C  1 

«23  '  12  NjhIEERS  FRC.T  UNIFCRtMO  7C  1) 

6TC  E I  -=  -  6 

650  FOR  I=,  TO  12 

o63  BI*e'Z»:<NO<0) 

670  NEXT  I 
6  79  ' 

675  •  NOI3Z=N<*1EAN»SE£U22) 

680  NCISZ-SO£',NEI*n£4N 

£35  ' - 

686  ' 

687  1  NCZSY  *IHZ  Dt.AY  MEASUREMENT 

6°C  7m£A3=70ELAY*NOISZ 

69  2  ' - 

£92  IFL  =  0  THEN  699  EuSt  700 

o9»  lPFInT  .!;R9X;9£Y;R&Z.RASRPSTi3£LAYiTnEAS 

7CC  - 

701  ' 

702  ' 

703  1  AUER  AGIO  TMtAS  AT  9. 5 . 29 . 5.  49. 5 , - .709.5  SEC. 

709 

710  S7hEAE  =  STM.:AS-h7rt£AS 
7.1  L*L-1 

718  '  CCH-L7E3  ANO  STORES  AUE RAGED  TKEA5  AND  RESETS 

u.STNEAS.ANj  INCREMENTS  N 

719  ' 

720  17  L*20  TrlEx  LxOlACIT'KM-STFEAS/ZO,  1N*N+1  ISTNEAS^j 
BE  1  NEXT  K 

37.  3 


HIS  PAG£  IS  Hi. 
JR<M  COil  FU.-.-l- 


3  lPR'INT'ESTIMATE  TARGET  PARAMETERS — MODIFIED  SAL  SS-NEF  TG.N  METHOD' 
i  caTaic .it.ic»j.9c .o.c.o.t 

1  3  ■ - 

EC  CIVEn  THE  ESTIMATE  Of  THE  TARGET  PARAMETERS.  TnlS  PRC  GRAM 
GENERATES  THE  JACEEIAv  AND  T(-£  FEE  CENG  FuNETIGn. 

3:  ■ - 

3:  ' 

33 

3i  '  TRANSLATION  OF  THE  NAMES  CF  THE  VAEIAEl.ES 

36  • - 

37  1  RS=SlANT  RANGE 

33  ' Ml =3uh IT  VAxEAEwE  USED  IN  LENE  1673 
A  D  '  AVTM= AVERAGED  TIME  DEhAY  MEASUREMENTS 
M  ' Ll=I7ERA TIv t  COjNTER 

IE  W1 ,*2.X  W3  STORE  INITIAL  X.Y.JZ  TARGET  POSITION 
H3  ' A-JAC03 IAN  KITH  -F<X)  IN  THE  SIXTH  COLUMN 
HR  '  C11  =  TAAMSFCR.MATZOn  ANG-E 
HI  P1I=FER-'CSMANU£  INDEX  EErGRE  HOUSEHOLDER 
H6  1 F2I  =  •  •  AFTER  ’ 

H7  U=USED  IN  HOUSE  HOLDER 
Ha  ‘ SIGN  =  •  • 

HO  'MAG=  *  ’ 

SC  SATA-  * 

51  ETA*  •  * 

IE  'X*  CORRECTION  TO  TARGET  PARAMETERS 
S3  6  =  USED  IN  SACK  SUE STETUTECN 

SH  'FI  AND  F2  ARE  FIRST  AND  L.AST  MEASUREMENT  WEIGHTS.  RESPECTIVELY 
SS  ' TT “TEMPORARY  STORAGE 

36  ‘ T  <  33 »  5 1  =  TEMPORARY  STO PAGE  FOR  INTERCHANGE  OF  RONS  CF  A 

57  ‘  KI= USED  IN  GliA  J  RATIO  FET 

58  '  ZI=USEO  IN  QUADRATIC  FIT 

59  '  IP=<KI>=  •  ■  •  • 

60  '  SA  =  •  ■  • 

61  '  AC*  •  •  .  .  . 

62  '  DETERMINATION  PARAMETER 

o3  '  02=  * 

6H  1  05=  OUTPUT  PE  TO  DISK  PARAMETER 
90  • 

no  ■ - 

1C  SO  ■  THE  FIRST  FIVE  CDl-MNS  CP  AOS. 5)  CONTAIN  THE  JACOBIAN 
WmElE  THE  SIXTH  <8.  LAST 3  COLUMN  CONTAINS  -F(X>. 

1010  ' 

1320  '  THE  AvTM ,  AVERAGED  TOE..AY  MEASUREhEHTS .  DATA  13  CALLED 
FROM  DISK. 

1C3C  DIM  AVTM(33).A<3S.3).U05).X(4).e'A).T<3i.5> 

1033  INPUT ’FILE SPEC 'I FILE SPEC* 
lO'-C  OPEN  .  1  .FILESFEC  s 
ICjO  FCF  L=0  TO  35 
106C  ENFcT^I.Av'TMIL) 

ic7c  print  ‘u  =  ' ;l; ‘avims1 ;avtm<l> ; 

1C  St  NEXT  U 
lO-'O  CLOSE 

1100  lPFENT* - 

lit:  LP  REN" 'FZLESP  EC  =  •  IFIuESF'ECR 

l 1 C  T  UP  PINT’ - 

11.3  • - 


VSlSPAS^ISBESTqfALin^-A 
;f  <  M  Iv.  I 


D-5 


X  *  < 


1 1 11  '  ***  3 :  T  TU *«  TC  WCIChT  T.-iE  FIRST  l  lAST  htASU.R  *  * £i-T  = 

i;io  i*f l r  * k r asu •.£»•! in r  wiiG.-iTi.vG'  >n-> 

ll3u  IF  Ni-  =  *NC*  Th EN  SI  TO  lloO 
l.*t:  INPUT • SET  FIR2  T  AnI  last  WEIGHTS  *  ♦  r  1 » r  2 
tlJJ  Lf  ftl.'  TBFI«ST  A.vl  LAST  WEIGHTS  ARC  J  •  JFl  #  'A/O*  5  FI 
n  6  j 

l ;  7,j  pri-t 

i  1 C  ?  FPIf*  a** **x  «**  ^ t»  «*»«.«  *»*»** «**:*«:«*  «/ixx  dixxx  «**;<,«  4C «#****«  • 

U^C  * 
nt:  • 

•  -  >  ■ _ , _ _ _ _ _ . _ _ _ 

1  IT  C  READS  I?  PUT  OAT  A  DEFINED  in  10CQ 

:i3.'  .-.saj  s:.s:»si.  jo/g:  »3r»it;L)f;;c.z2r 

l..*t  •*•  317  IM'l^  TARGET  F  AF Ar ITER  :  *  ** 


11*3 

1Z7C  TT  =  - 1 

nu: 

:z'-c  vxr-c:vY ~*c 

12  31  R'3=::cc  c 

131  i  R.Z7*3Q0C 

1310  *A=RS>  *40 0 1  3  ) 

III }  Rr-TT*(KS-<.':S-FA«KA-PZT*r.ZT)C  .5 

13*1*  C11  =  <G0-1S.  )  f  .tf:?4E22 

1 3i c  :*cc3<c:1)*fa«s:'  <c:n 

13.il  ftYTa-RTaSlMCil  >  *.<A*C  J3  (  C 1 1  ) 

12aj  31=0 
2  37C  ’ 

132:  1  ***  LL=ITEF-ATI^'E  CDJNT  **# 

1 3  7  3  ' - 

1ACC  Ll-lL+I 

1-01 

i-i :  Ki»o:zi*a 

1413  1*1 

1*1?  IF  r.I=  ?  TnE.v  FF:lNT,WrVCNG  DIPS  I TICN * J  £Nl* 

141=  • 

;-2L  lf  c  »•«**-*»  mh< 

144;  ^;PINT*C6SL  AvIR  AnI-  TARGET  PARAPETS  \S  *  *  'Ll-  *  *  Ll  #  *KI=  *  »KI »  *  Z'<  ( Kl-i 

I-: : ;  ■zi-**  ;zis  fSA«*  ;sa 

14  5C  _ PRINT  * - • 

14cl  w.  r.lNTS  j*SltS2. 33»Gl  >  G2  »  ZCC  ,  Zl  3 »  Z2G ♦  </XT  f  VY  T » R/T  ,  ftY’ » RZT 

l*7l  SETf  CtVIRvEF:  INITIAL  SPIED  AND  HEADING  IN  RADIANS 

1  4  o  ;  _r  .\  !?■  T 

14-::  w  f  r.lN*  ' - • 

1  SC  i-  1  ***  ST  Or  ..a  TARGET  RANtJI  CCH -OnENTS  «*< 

is;-  wi  —  a  :*2  =  <y~  :w  3^r~ 

if.  It  20 -S*J  :GAftlC*££  *.  01  7*53'J 

1S2?  vlXJ-SC^SlN  GAlAO)  :  vIvC*=('4iCO:  '  GAMA!)  JUZO-ZCC 

lS^  I'ITIAw  F  I .» I T  I L  N  Cr  T.-tl  Cl  SliVlr  #I.E.  *FC3ITIN  AT  ?.S  SECONDS 

:  l  i  c 

1Sol  ? 1  XT  =  *  1  a  .  «V  .  3  5RLY3=vZVO«?*5  !k32C«VZ  3«?  .  3 

;:?C  ’Ill "ION  Cr  The  TAhCETiaT  7 . 3  SECONDS' 

1339  *«*  STD  - 1  IhlTIA:.  TARGET  F35I"I3n  **« 

1  S^C  F>  :  -y<T**.  ->p  XT  :  *,  r  '  -  ^  >  T  « ’■* .  S  Y  T 

louO  ' - * - - - - - - - 

lc>t  RELATIVE  PCIITIln  AT  9.S  SLIIUS 

1  6 1  l  ■  4> -X  'T-R  L  *  .  :s<  jY  -r»  T  T-F  I YI  J <“.!  =RZT  -k3Zj 

1  L  3  j  fv- -'Fi/f'-  •  *  .  f  •  i  :■  r  •  ivi  I  •  F  ol ) .!  .  3 


i  S  PASt 


■:.i  ;■: 

1  j . . 


D-6 


.  x 


1 

! 


i 

i 

? 

1643  • 

1 65  C  *  — - — — — — - - -  —  —  —  —  — - - - . m - _  _  _ 

16o0  '  -F  <  0  )  AND  u(0.1,2.2.9>  AT  TIME=9.5  SECONDS  STORIED 

IN  A(O.I).  F  ANO  J  DON'T  APPEAR  IN  THE  PRCCRAM.  *xx 
1670  Cl l=GAMAO-ltfl  «  .  0179322  :RA«n9:<*SIn< Cl  1  )-fR5Y*C3S< Cl  1  ) 

1680  ' 

1690  A<0 .3) =■- (RA-RSx^OOORAvTMC  0 >  J 

1703  A(3<9)«-9CaC»AvTM<C)  *.<4Z/RS :  A<  0  «3>*CCS  (211)-9QOOxAVTM<G  )»R5Y/PS: A<  C  »2)*SIN< 

Cil )~90J0»AVTM<O>»R9X/RS 

1710  A(0f  l)-9.5xA<0»3):A(0i0)  =  <:.S*A(0i2> 


1793  FlI»A<a.5)*A<0.5) 


:  /  j.  ■ 

1760  '  MAIN  LOOP  K*1 <  TIME*29 . 5)  TO  K»33  <TIMt=/C9.5  SEC.) 

1770  1  »»************x*»***:«**x****«***x******  jiuiiiimii 

1783  • 

1700  ■ 

18C3  FCF:  K=1  TC  33 

1810  '  CHECKS  FOR  CESERvER  MANEUVERS 

1820  IF  K=»12  THEN  1830  ELSE  1950 

1830  •  COMPUTE  POSIT  ION  AT  22°  SECONDS .  THEN  LF  OAT E  SO  S  GAMAO 
1890  ' 

1833  R1XO«V1XO»9.3+R'1XO  SR2rC=VZY0*9 . 3-R2YC  :(}?ZC=VZ099 . 5  +  R3ZC 
I860  > 

187C  S0=S1  :GAMA0-G1». 0179333  1020=210 
1880  VIXO-SOxSInCGAmAO)  1 V2YO=SQ«CCS< SAmaC > 

1E90  ■ 

1 °  3  0  '  **»C3MPUTZ  C6 SERVER  POSITION  AT  K*12<  TIME *2^9 . 5  SEC.)  »** 

19  10  R1XC-'V1XC*1 0 .5*R1XO;R2YO»</ZYO*13.5+R2Y3.M2ZO-VZO*1  3 . 5+R32J1G0TQ  2110 

;92Q  ■ - 

1 93C  ' 

179)  '  CHECKS  FCr  SECOND  OE SERVER  MANEUVER 

1930  IF  K=Z9  THEN  1973  E_SE  2070 

19*0  '  COMPUTE  QB3.  PCS.  AT  979  SEC.  THEN  UPDATE  SO  t  GAMAC 

1970  R1X0=V1XC«9.5+R1X3  1R2YC=V2Y3»9 . 5>R2YO  1R3ZC=VZ0»9 . 3<-R2ZC 
1980  ' 

1900  SC=S2  1GAMA0=GZ*. 0179333  1VZ3*Z23 
200  0  V1XC»SC*SIN'.GAMA0)  .‘VZYO=SO»COS<GArtAJ> 

20.0  '  COMPUTE  033.  POS.  AT  K=29 < TIME=989 . 5  SEC.) 

20  20  ' 

2C30  R1XO=V1XO»10  ,S+RlX0t?2Y0*V2''0*lC  .5»R2Y0:R3ZC*VZ0*ia  .S+R3Z0 

2093  GCTO  2110 

2C33  1  IF  A  MANEUVER  CCCLR5  THE  NEXT  POSITION  IS  COMPUTED 

AECVE  .  SC  THE  LCCIC  ShCvN  EEuCH  IS  SKIPPED. 

2060  1 - 

I'./O  1  X.Y,  AND  2  POSITION  07  THE  OBSERVER 
2  3  E  3  ■ 

239)  RiXO=VlxQ*20*R:XO  !R'2YO-v2YC*23*ft2YC  lR2ZO=vZI*2Q*ft3ZC 

2110  1 

2121  'ESTIMATE  OF  X,Y.  AM  2  POSITION  CF  THE  TARGET 

2130  Rz;»UYT*209kXT  JRYT=VYT*20+RyT:R:ZT-R'2T 

2 :  9C  •  *»*  ES'lMA'E  OF  RE.ATIVE  FCEI'ICN  •««* 

21  jf  ’ 

2I6C  F  9»-kXT-R  ixc:  R3Y=PYT-R2Y0:  RSI  -RZT-R30 
2 . ’3  R 3- i R6X»R9>  *R3Y  »R5f ♦AoI«  <42)0 .5 

2131  ' - 

2;°o  •  «**  elen.m's  3F  the  f-array  stc-ed  as  -f  in  aok.si  <**  *** 


D-7 


•  i  xci 


~  ^4 J 


zrac  • 

2210  Cll«:ArA:-131«.Cl^'3: 

222  J  RA-StX«SIn*  C 11 
223C  rtl*403?*A  /^a') 

22*1  A(K»  5;  *- <P.A*-n£*h l  ) 

223  c  ‘ 

226C  '  *x*  £.£hiiNTS  C.“  THE  J-AfcFAr'  *** 

• _ _ _  ^ _ 

2269  A<K,9>=-I».* R6Z,  P.S 
2293  Air  *  3  )  =20i  '  C.  1  >  -Hi*  12T/.3S 
25  J  j  A  ( .1 » 2 ) *SZN i C 1 1  >-M:«R*t.t/AS 
"313  =  '  «A<;  »3< 

2333  Avl  fC)  =  (;I*2I  — ‘F*w.'*Af.  K  r  2 ) 

323'.  ' 

32H  F-1I=F‘1I-A(. 3) 


23ot.  ' - 

236-  '  ***  OLTFuT  fT  TO  32S.I  *»» 

2269  2E  03  =  1  T.  CM  2 COO 


236o  If  U.=  0  THEM  2323 

23-63  ’  *<«  T2,l.n:;-AT  I  EV  INECT.—ING  true:  target  F  ARA.iXTGRS  «** 

22*i  G1  =  (.:-1I-F3I'/F23  :OZ=AE:StC  1> 

2369  If  Q2  =’,1  l»tN  93/A 

23/9  ' 

23/3  1  *»«  S-A9T  Cr  GUA..DATI  2  FIT  *** 

2332  IFF <K1>»F II 

2333  IF  LL=0  THEN  22": 

23"C  IF  ZI=2  THEN  23/0 
2*-  J G  IF  ZZ-\  THEN  253C 

2910  IF  IF-r'DI;  =F’I'i_i.-n  THEN  2“50 

2920  SA-— 1/(21  CD 

293E  322L0;  M5C 

2/90  GCDC  192C 

295'.  IF  ill  l  TliC.9  2990 

2960  SA=i:ZI=l 

2970  GC3UJ:  9-50 

2a;o  GQrC  1920 

299  0  A2=<3».-:<Li_-l  )-9»I--FC:I)*IFF'<KI-l  )  )/(2*FI<Ll.-l)-9«Ir-r-<KI  ■♦2»I"F  ( *1-  1  )  >  !  ZI  =  2 
2500  SA  =  <  A2-l)/(2:  ull-l  >  J 
23  10  GOEL'J  9*3  1 
2320  GET 3  .920 

2’3C  IF  IFF'  ( ,:I )  IFF  (,:2-l)  Then  25/: 

2291  EA  =  .  -F  1  -9.1;  MM  )-3«I;  iHI).,  <  2*(  F  I  (  LG-  1 ) -2  'IP  -  ( FI- 1 '  *  IF  F  >  ,;i )  )  1 1 Z2  r 

22  j  i  OC  ' L  t  4  DC 
226C  GOT"  i*r: 

25" 0  ri<LL)*:f-t  'CD 

2273  '  »**  ENG  O.-  3. Ah.  ADO  FIT  »*« 

2332  lF  •  2 n T  • - - - - - ■ 

*.  2  3l  Lr  I-  2*rT  ‘  !  .1  i  *  r  1  I  Ol  -  *  I  (  I  (  1 _ J 

ZED  :.:.3  =  6  :  D--.-3.-; 

25"- :  :r  u*  >  t.-d  t  • 

IF  L-«9  ’  3D. 

2 i 0  o  -  - - - - — - - - - - 

2ajg  «■»  hi_:o,,("  ir.j  g  ri:  :i  s  ^2  rino-j >  :/Et  ;  »»» 


( 

l 

I 

J 


] 

1 


0-8 


’IS  PACT  IS  BEST  QUALITY  PRAGUE** 

•’  uti-x  i\wi.  1...L  ;  a  udc 


2413 

2oZ0 

2630 

264C 

265C 

2660 

2664 

2665 

2670 
2a£C 
2690 
27  3  3 


2724 

t  ~  ^ 

272  J 
2B4t 
2B3C 
2860 
2870 

28  JO 

2ar-o 

29  C  0 
29  ;C 
2°20 
2930 
29-H 
2°50 
296  0 
2970 
2983 
2490 
20  C  C 
3C  C 
3070 
3380 
3C  °0 
3100 
3120 
314C 
3130 
3 1  6C 
3170 
3LuO 
31  40 
22  0  0 
31 1 C 
3220 
3220 
324  0 
3300 
3310 


IF  N»='N3'2MEN  GC'O  2710 

lFRINT'AI  C  .5)  =  '  iA'.  0.1)  !  'A<25.5)*'  ;a<33.3> 
FOR  1*0  TC  5 

A<0.I)=F1«A<0,I) 

A<33.I>»F2*A<35.I> 

NEXT  I 


•  in  3ECAUSE  THE  FIRST  t  lAS”  ROWS  HALE  THE  LARGEST  E-EMENTS 
WHEN  WEIGHTED.  TH£>  ARE  STACKED  TOGETHER  **» 

'  INTE  :.CHAnG£  THE  OTH  A  M2  THE  34~H  ROWE  OF  A(K.J) 

FOR  U*C  TO  5 

TT=A<34.d) :A<34.J)*A(0.J) :A(0.J)«TT 
NEXT  J 


'  in  BECAUSE  THE  LAST  ROWS  07  MATRIX  A  HALE  THE  LARGEST  ELEMENTS 

THEY  ARE  S 'ACTED  FIRST  FOR  THE  HOJSEHOlDER  TRANSFORMATION  in 
1  INTERCHANGE  THE  ROWS  CF  A 
SOS JE  3380 

'  ***  CAlL  householder  *** 

GO SUE  3230 


ER!NT't«***»***M*««*9*t«***««M**«***t«««4»**t*»*#M*****4*' 


'»»«  BACK  SUESTITJTION  IS  NEXT  *»« 

FOR  1*3  TO  4 

e<I)*A(I,5) 

NEXT  I 

FOR  L*4  TO  0  STEF  -1 

X(L ) *E < L ) /A( L  L ) 

FOR  1*0  TO  L-l 

B(I)-E-<I)-A<I.L)«X(l) 

NEXT  I 

NEXT  L 


'  ***  UPDATE  TARGET  PARAMETERS  *»« 

WXT *UXT*X<  3  >  ;VYT*DYT*X( 1 ) 

RXT*witX(2)  :RYT*h2*X<2)  :RZ',*w3*X(4) 

LRRIHT' - 

LF  RINT 'LL*  *  ILL i *  X<3)  THRU  X(4)  • SX< 0 ) i X( 1 ) i X( 2 ) t X<3> SX ( 4 > 

LFRIfll  * - 

'  93(«43&Cni  END  OF  MAIN  LCCF  990l?9'»O9ge 

GOTO  1400 


***  SUE RC jTIf*£  HOLSEHCwDER  *»« 

**«  mOUEEmSI-DEF  TRANSFORMATION  *«* 

RED JIRED  INFL*  .CELS  SRO»St  CCLS-f  RINT  .  ROWS*  "RINT  ) 
NET  REQuI- ID  IF  CON1AINED  IN  MAIN  PROGRAM 


0 


»n*ClL5-  1 


T 


•  JO 
2  ;  - 


SET  T!iZ  NL'f  L  £  f\  Gr  PJ-  ;  .  C «_  L  3 
l  M  -  *  l*  A  5  -  1 

**»  rtAIN  ^O'Jr  t  TRAv::->ttS  A  A ~  A 

STO  '■£  ft’  CCu- <  N- 1  ) — AG 3  jm£3  LAST  C C  . 
FOA'IlMi  rC.'<C  “ION.  *** 


ri.‘.L . 

IE 


fc  *»  k«: 


349 : 

24 2 

3«t2J 

2 - i-9 

3- :o 


ora* 

3— : 

3iC  J 

:zzi 
352C 
254  *J 


2'lZO 
35*  £ 
32  r: 
2523 
25  >9 
:* ): 
36 :  J 

2t  *. : 


3**0 
3oE ./ 

Jo  r  v 


AFP  A  i  <  ** 


:c  ( N- l / 

MA’^*0 

**«  CW-TCS  HAUMTlOE  CP  55-  r 
FC*  :*K  TC  M 

IF  K*h  THEN  30 TC  3Bo3 

***  the  Atc-.t  c.-»£r::s  p: 
MAG-M  43-h  (.  I  f>'- )  *A  K  Z  ) 

NEXT  Z 

***  SE’3  N£V*  DIAGONAL  ELEMENT 
SIGN* A  ( 1. 1 K  )  /ACS <  A  ( » K  )  ) 

MAG*-SIGN«'MAGC .5) 


CONFUTES  U  YEITCP 
U(K)«A<K* 'O-ttAG 
FO.\  w«K*i  TC  M 

Ui J)aA( J»K  ) 

NeXT  j 

MMX  K-TH  CO_  13  LrOATID  STARTING 
AT  THE  K-TH  F:3*.  RCP  ELEMENTS 
FRC.1  K  +  l  On  APE  SET  TC  0(SEE  2180) 

A(K»K>*i1AC- 

XXX  THE  FENhINInG  K-l  CClS  APE  NOW  T = ANiF JRPEO 

BA*A»l/'MA3*-..0:>  > 

FC*  w-=K  +  :  TC  N 
ETA-C 

FC*  I-K  TO  n 

ETA=ETA  +  'J<I>*A(I,w> 

next  I 

eta*-:ta*c.ata 

***  ELEMENTS  of  the  I-TH  CCl  OP  a  AkI  PZt'EF :>*£:•* 


•  A '  J  *  K  >  =  C 


FC*:  I-K  TO  M 

s*: ;  a  ( r » tj  >  =  a<  i  •  w  )  +e*a  <  j  ^  i  * 

271-  NE  <  r  I 

3'C:  **«  rtJXl/v?;  TC  THc  NOT  COuJiN  *•<« 

37  5c  NEXT  sJ 

3  7  ■•  C 

37::  moving  rn  the  not  Tkansfoematio^  *** 

C7TI  Nt  T  K 

S7o  0  F:‘P.  I  ..'.•£■/  T'4C  -  SO.EF.-  AFTZP  Tf*  An  ;r  £.FM.-»7IC?  *** 

37  ; •*  f  II-PCT'-A ■  !•:  > 

3iw)  Nt  T  : 

r  5 -f  r  IN  T  1  *  <11  ;«!<0<;4t  »t*f  * 


*  w“* 

...  l»  4A*J* 


D- 10 


-i 

J 

J 

1 

] 

1 


“G3:  • 

icnl  "LL*#  #Ll. 

3B30  1 

33Au  LPRIN7 *«M»# <4 #*«♦#♦*#** 

:er:  return 

387.! 

31*76  ’ 

333C  1  SLSRCjTXNE  PC  A'  INTE  v  JmAnGInG  ALL  ROWS  OF  A(35>5:.  RETLT  MS  TI  2353 

3 8>5  FDR  K«0  TC  33 

2<-::  F C.<  J*C  TC  5 

39. C  T(K*J)»ACKf J) 

3°C0  NEXT  J 

3930  NEXT  K 

39 4 Q  FCR  :-(*0  TC  33 

3950  FOR  j»0  TO  5 

3960  A(33-Kt  J>*T(,<fJ) 

3Q7 0  NEXT  o 

396?  NEXT  K 
3990  RETURN 

*J0C  - 

4570  - - 

4273  ’  ***  FROM  A  GOTO  IN  LINE  2369  *** 

4374  PXUL>»PlI 

4375  LF?ttNT*U.*“  5 LL • " F I (Lu ) - '  JP1X. "GZ-* iQ2* •PZX-*  tF2r 
4373  Q5-1 

42UC  IM>uT#,JXT*#  5LXT 
4391  INPUT  *  VYT  =  *  #  VYT 
44C3  INF  JT  ,RXT**JhXT 
4410  INPUT  #RYT*'5PYT 
4^20  INF  UT  *R2Ta • » RZT 
4430  GOTO  1400 

4445  1 - - - - - - - — * - * - 

**,50  '  «**  SUilOuTIvE  TO  UrCAT!  TAf-SET  PARAMETERS  »■** 

4Ao0  '  CALLED  <  lINEE  2A3a.2Ar0.25ie.«2250  >  EY  QUADRATIC  FI* 

<M*C  vXT-VXT*SAxX(U)  !VYT=viYT*SA*Xll)  !  RXT*m  +SA*:<  ( 2 )  tRYT*K2*SA*<(3>  !RZT=w3*Sa 
*X  (A) 

AA8C  RETURN 

‘('IRC  - - 

5300  '  OUTPUT  TO  OI3r( 

50 C  l  '  *»*  FRCM  A  5CTO  IN  LINZ  236A  «** 

5305  PKLU-F’lI 

5  3  C6  LPPINT  •L-»'JLL.*PI<L->«*.‘PI<LL> 

5313  Of£N,O,»2»,PI<LL)/0U:i* 

5013  FOR  • _ =  3  TO  9 

503t  PRINT *2  >P  I( -L  i 
5C  AO  NEXT  LL 
5050  CLOSE 
SieO  ENO 


0-11/0-12 
Reverse  Blank 


INITIAL  DISTRIBUTION  LIST 


Addressee  No. 

NAVSEA  (SEA-63,  -63R,  -63R-1,  -63R-13,  -63X,  -631X 
-63D  (2  cy) ,  PMS-409  (2  cy ) ) 


OTIC,  Alexandria 


of  Cop 
10 

12 


END 

DATE 

FILMED 


I 

* 


DTIC 


